from __future__ import division

################################################################################
#
# convert_date_time.py - Provides date/time conversions as part of the 'convert'
#                        subpackage. Not to be imported by itself!
#
################################################################################

__author__ = 'Adrian Price-Whelan <adrn@nyu.edu>'

# Standard library dependencies
import re, sys
import math as m
import datetime

# Third-party
import numpy as np

# astrony
import geometry

def parseRA(raString):
    """ Parse an RA string and return hours, minutes, seconds.
    
    Parameters
    ----------
    raString : string
        A string representing a Right Ascension (RA). The values can be delimited by
        white space, no space, or a colon:
        13:19:23.1 or
        13 19 23.1 or
        131923.1
    
    Returns
    -------        
        
    """
    print raString
    div = '(|:|\s|/)'
    pattr = '([+-]*\d?\d)' + div + '(\d?\d)' + div + '(\d?\d[\.0-9]*)'
    try:
        elems = re.search(pattr, raString).groups()
    except:
        print "Invalid input string!"
        raise

    hours = int(elems[0])
    minutes = int(elems[2])
    seconds = float(elems[4])
    
    # Check for nonsense values
    if hours > 24.0:
        raise ValueError("Hour value must be < 24.")
    if minutes >= 60.0:
        raise ValueError("Minute value must be < 60.")
    if seconds >= 60.0:
        raise ValueError("Second value must be < 60.")
        
    return (int(hours), int(minutes), seconds)
    
def parseDec(decString):
    """ Parse a Dec string and return degrees, minutes, seconds.
    
    Parameters
    ----------
    raString : string
        A string representing a Declination (dec). The values can be delimited by
        white space, no space, or a colon:
        +13:19:23.1 or
        -13 19 23.1 or
        +131923.1
    
    Returns
    -------        
        
    """
    
    div = '(|:|\s|/)'
    pattr = '([+-]*\d?\d)' + div + '(\d?\d)' + div + '(\d?\d[\.0-9]*)'
    try:
        elems = re.search(pattr, decString).groups()
    except:
        print "Invalid input string!"
        raise

    degrees = int(elems[0])
    arcminutes = int(elems[2])
    arcseconds = float(elems[4])
    
    # Check for nonsense values
    if degrees > 90.0:
        raise ValueError("Degree value must be <= 90.")
    if arcminutes >= 60.0:
        raise ValueError("Arcminute value must be < 60.")
    if arcseconds >= 60.0:
        raise ValueError("Arcsecond value must be < 60.")
        
    return (int(degrees), int(arcminutes), arcseconds)

# BELOW HERE NEEDS WORK!
def raLST2ha(ra, lst, raUnits='hours', lstUnits='hours', haUnits='hours'):
    """ 
    Converts a Right Ascension to an Hour Angle, given the Local Sidereal Time.
    
    Parameters
    ----------
    ra : int, float, Angle, string
        Accepts a Right Ascencion as any of the above types.
    lst : int, float, Angle, string, datetime
        Accepts a Local Sidereal Time as any of the above types.
    raUnits : string
        Can be either HOURS or DEGREES, defaults to HOURS.
    
    Returns
    -------
    hour angle : Angle
        Returns the computed Hour Angle as a geometry.Angle object
    
    Examples
    --------
    >>> raLST2ha(70.3425, hour=14, minute=26, second=18, raUnits='DEGREES')
    9.74883333333
    >>> raLST2ha(5.055477, hour=14.4383333333333333)
    9.74883333333
    
    """
    raAngle = geometry.Angle.fromUnits(ra, units=raUnits)
    
    if isinstance(lst, datetime.datetime.now().__class__):
        lstAngle = geometry.Angle.fromDatetime(lst)
    else:
        lstAngle = geometry.Angle.fromUnits(lst, units=lstUnits)
    
    return lstAngle - raAngle

def haLST2ra():
    pass

def raUTC2ha():
    pass

def haUTC2ra():
    pass

def ra2transitTime():
    pass
    
def raDec2AltAz(ra, dec, latitude, longitude, datetimeObj, longitudeDirection='W'):
    """ 
    Converts an RA and Dec to Alt Az for a given datetime object.
    
    Parameters
    ----------
    ra : float (any numeric type)
        A Right Ascension, default units DEGREES
    dec : float (any numeric type)
        A Declination, default units DEGREES
    latitude : float (any numeric type)
        A latitude in DEGREES
    longitide : float (any numeric type)
        A longitude in DEGREES
    datetimeObj : datetime.datetime
        A Python datetime.datetime object
    longitudeDirection : string
        Default is longitude WEST, 'W', but you can specify EAST by passing 'E'.
        
    Returns
    -------
    alt : float
        Altitude, default units DEGREES
    az : float
        Azimuth, default units DEGREES
    
    """ 
    gmst = utcDatetime2gmst(datetimeObj)
    h,m,s = gmst2lst(longitude, gmst.hour, gmst.minute, gmst.second)
    lst = (h + m/60.0 + s/3600.0) * 15.0
    ha = lst - ra
    
    alt_rads = asin(sin(radians(dec))*sin(radians(latitude)) + \
        cos(radians(dec))*cos(radians(latitude))*cos(radians(ha)))
    
    cos_az = (sin(radians(dec)) - sin(alt_rads)*sin(radians(latitude))) / \
        (cos(alt_rads) * cos(radians(latitude)))
        
    alt = degrees(alt_rads)
    az = degrees(acos(cos_az))
        
    if sin(radians(ha)) < 0:
        return alt, az
    else:
        return alt, (360.0 - az)
    
def eclipticLatLon2RADec(lat, lon, latLonUnits='DEGREES', raDecUnits='DEGREES'):
    """ 
    Converts an Ecliptic Latitude and Ecliptic Longitude to a Right Ascension and 
    Declination.
    
    Parameters
    ----------
    lat : float (any numeric type)
        An ecliptic latitude, default units DEGREES
    long : float (any numeric type)
        An ecliptic longitude, default units DEGREES    
    latLonUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
    raDecUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
        
    Returns
    -------
    ra : float
        A right ascension, default units DEGREES
    dec : float
        A declination, default units DEGREES
        
    Examples
    --------
    >>> eclipticLatLon2RADec(70.3425, -11.4552)
    (143.72252629028003, 19.535734683739964)
    
    """ 
    if latLonUnits.upper() == 'HOURS':
        lambdaa = lat * 15.0
        beta = lon * 15.0
    elif latLonUnits.upper() == 'DEGREES':
        lambdaa = lat
        beta = lon
    else:
        raise AssertionError('latLonUnits must be either HOURS or DEGREES')
    
    # Python works in Radians...
    lambdaaRad = lambdaa * math.pi / 180.0
    betaRad = beta * math.pi / 180.0
    eeRad = ee * math.pi / 180.0
    
    deltaRad = math.asin(math.sin(betaRad)*math.cos(eeRad) + math.cos(betaRad)*math.sin(eeRad)*math.sin(lambdaaRad))
    
    y = math.sin(lambdaaRad)*math.cos(eeRad) - math.tan(betaRad)*math.sin(eeRad)
    x = math.cos(lambdaaRad)
    
    alphaPrime = math.atan(y/x) * 180.0 / math.pi
    
    if y > 0 and x > 0:
        alpha = alphaPrime % 90.0
    elif y > 0 and x < 0:
        alpha = (alphaPrime % 90.0) + 90.0
    elif y < 0 and x < 0:
        alpha = (alphaPrime % 90.0) + 180.0
    elif y < 0 and x > 0:
        alpha = (alphaPrime % 90.0) + 270.0
    else:
        alpha = alphaPrime
        
    if raDecUnits.upper() == 'HOURS':
        ra = alpha / 15.0
        dec = deltaRad * 180.0 / math.pi / 15.0
    elif raDecUnits.upper() == 'DEGREES':
        ra = alpha
        dec = deltaRad * 180.0 / math.pi
    else:
        raise AssertionError('raDecUnits must be either HOURS or DEGREES')
    
    return ra, dec

def raDec2Galactic(ra, dec, latLonUnits='DEGREES', raDecUnits='DEGREES'):
    """ 
    Converts a Right Ascension and Declination to an Galactic Latitude 
    and Longitude
    
    Parameters
    ----------
    ra : float (any numeric type)
        A right ascension, default units DEGREES
    dec : float (any numeric type)
        A declination, default units DEGREES    
    latLonUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
    raDecUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
        
    Returns
    -------
    b : float
        A Galactic latitude, default units DEGREES
    l : float
        An Galactic longitude, default units DEGREES
    
    """    
    if isinstance(ra, [].__class__) or isinstance(dec, [].__class__):
        ra = np.array(float(ra))
        dec = np.array(float(dec))
    elif isinstance(ra, np.array([]).__class__):
        ra = ra.astype(float)
        dec = dec.astype(float)
    else:
        ra = np.array([float(ra)])
        dec = np.array([float(dec)])
    
    if raDecUnits.lower() == 'degrees':
        ra *= toRad
        dec *= toRad
    
    raNGP = 192.85948 * toRad
    decNGP = 27.12825 * toRad
    lASCEND = 32.93192 * toRad
    
    gb = np.arcsin( np.cos(dec)*np.cos(decNGP)*np.cos(ra - raNGP) + np.sin(dec)*np.sin(decNGP) )
    gl = np.arctan2( np.sin(dec)*np.cos(decNGP) - np.cos(dec)*np.cos(ra - raNGP)*np.sin(decNGP),\
                     np.cos(dec)*np.sin(ra - raNGP) ) + lASCEND
    
    gl /= toRad
    gb /= toRad
    
    gl[gl < 0.] += 360.
    
    if len(gl) == 1:
        return gl[0], gb[0]
    
    return gl, gb
    
def galactic2RaDec(gl, gb, latLonUnits='DEGREES', raDecUnits='DEGREES'):
    """ 
    Converts a Galactic Latitude and Longitude to Right Ascension and Declination
    
    Parameters
    ----------
    gl : float (any numeric type)
        An Galactic longitude, default units DEGREES
    gb : float (any numeric type)
        A Galactic latitude, default units DEGREES
    latLonUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
    raDecUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
        
    Returns
    -------
    ra : float
        A right ascension, default units DEGREES
    dec : float
        A declination, default units DEGREES
    
    """
    
    if isinstance(gl, [].__class__) or isinstance(gb, [].__class__):
        gl = np.array(float(gl))
        gb = np.array(float(gb))
    elif isinstance(gl, np.array([]).__class__):
        gl = gl.astype(float)
        gb = gb.astype(float)
    else:
        gl = np.array([float(gl)])
        gb = np.array([float(gb)])
    
    if latLonUnits.lower() == 'degrees':
        gl *= toRad
        gb *= toRad
    
    raNGP = 192.85948 * toRad
    decNGP = 27.12825 * toRad
    lASCEND = 32.93192 * toRad
    
    dec = np.arcsin( np.cos(gb)*np.cos(decNGP)*np.sin(gl - lASCEND) + np.sin(gb)*np.sin(decNGP) )
    ra = np.arctan2( np.cos(gb)*np.cos(gl - lASCEND), \
                     np.sin(gb)*np.cos(decNGP) - np.cos(gb)*np.sin(gl - lASCEND)*np.sin(decNGP) ) + raNGP
    
    ra /= toRad
    dec /= toRad
    
    ra[ra > 360.] -= 360.
    
    if len(ra) == 1:
        return ra[0], dec[0]
    
    return ra, dec

def raDec2EclipticLatLon(ra, dec, latLonUnits='DEGREES', raDecUnits='DEGREES'):
    """ 
    Converts a Right Ascension and Declination to an Ecliptic Latitude and 
    Ecliptic Longitude.
    
    Parameters
    ----------
    ra : float (any numeric type)
        A right ascension, default units DEGREES
    dec : float (any numeric type)
        A declination, default units DEGREES    
    latLonUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
    raDecUnits : string
        Can be either HOURS or DEGREES, defaults to DEGREES.
        
    Returns
    -------
    lat : float
        An ecliptic latitude, default units DEGREES
    lon : float
        An ecliptic longitude, default units DEGREES
    
    
    """
    if raDecUnits.upper() == 'HOURS':
        ra = ra * 15.0
        dec = dec * 15.0
    elif raDecUnits.upper() == 'DEGREES':
        ra = ra
        dec = dec
    else:
        raise AssertionError('raDecUnits must be either HOURS or DEGREES')
    
    # Python works in Radians...
    raRad = ra * math.pi / 180.0
    decRad = dec * math.pi / 180.0
    eeRad = ee * math.pi / 180.0
    
    betaRad = math.asin(math.sin(decRad)*math.cos(eeRad) - math.cos(decRad)*math.sin(eeRad)*math.sin(raRad))
    
    y = math.sin(raRad)*math.cos(eeRad) + math.tan(decRad)*math.sin(eeRad)
    x = math.cos(raRad)
    
    lambdaPrime = math.atan(y/x) * 180.0 / math.pi
    
    if y > 0 and x > 0:
        lambdaa = lambdaPrime % 90.0
    elif y > 0 and x < 0:
        lambdaa = (lambdaPrime % 90.0) + 90.0
    elif y < 0 and x < 0:
        lambdaa = (lambdaPrime % 90.0) + 180.0
    elif y < 0 and x > 0:
        lambdaa = (lambdaPrime % 90.0) + 270.0
    else:
        lambdaa = lambdaPrime
        
    if latLonUnits.upper() == 'HOURS':
        lat = lambdaa / 15.0
        lon = betaRad * 180.0 / math.pi * 15.0
    elif latLonUnits.upper() == 'DEGREES':
        lat = lambdaa
        lon = betaRad * 180.0 / math.pi
    else:
        raise AssertionError('latLonUnits must be either HOURS or DEGREES')
    
    return lat, lon
    

if __name__ == '__main__':
    import unittest
    import numpy as np
    
    class TestConvertCoordinates(unittest.TestCase):
        def test_parseRA(self):
            for h,m,s in zip(np.random.random(1000),np.random.random(1000),np.random.random(1000)):
                self.assertEqual(parseRA("%02d:%02d:%05.2f" % (int(h*24.),int(m*59.3),s*59.3)), parseRA("%02d %02d %05.2f" % (int(h*24.),int(m*59.3),s*59.3)))
                self.assertEqual(parseRA("%02d:%02d:%05.2f" % (int(h*24.),int(m*59.3),s*59.3)), parseRA("%02d%02d%05.2f" % (int(h*24.),int(m*59.3),s*59.3)))
        
        def test_parseDec(self):
            for d,m,s in zip(np.random.random(1000),np.random.random(1000),np.random.random(1000)):
                if (d+m+s) > 1.5: pm = "+"
                else: pm = "-"
                self.assertEqual(parseDec("%s%d:%02d:%05.2f" % (pm, int(d*90.),int(m*59.3),s*59.3)), parseDec("%s%d %02d %05.2f" % (pm, int(d*90.),int(m*59.3),s*59.3)))
                self.assertEqual(parseDec("%s%02d:%02d:%05.2f" % (pm, int(d*90.),int(m*59.3),s*59.3)), parseDec("%s%d%02d%05.2f" % (pm, int(d*90.),int(m*59.3),s*59.3)))
        
    unittest.main()


